RMSE & Bias
cols <- c("dataset_id", "formula", "shape", "fit_link_transformed")
df <- result_df %>%
filter(fit_link_transformed == data_link_transformed) %>%
mutate(
fit_link_transformed = factor(fit_link_transformed, levels =
c("Log",
"Softplus")),
data_link_transformed = factor(data_link_transformed, levels =
c("Log",
"Softplus"))
)
# full RMSE results
p_rmse <- df %>%
ggplot(aes(y = fit_family_transformed, x = rmse_s)) +
coord_cartesian(xlim = c(0, 1)) +
scale_x_continuous(labels = c("0", "0.5", "1"), breaks = c(0, 0.5, 1)) +
geom_boxplot() +
stat_summary(fun="mean", shape = 4, size = 0.3) +
facet_grid(data_link_transformed + shape ~ data_family_transformed) +
xlab(TeX("$RMSE(beta_{xy})$")) +
ylab("Fit likelihood") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
axis.ticks.y = element_blank())
p_rmse
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).

ggsave(paste0(figure_path, "positive_rmse.pdf"), width = 210, height = (297/4) * 3, units = "mm", useDingbats = TRUE)
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
# RMSE log showcase
scales <- list(
scale_x_continuous(limits = c(0, 0.2), breaks = c(0, 0.1, 0.2), labels = c("0", "0.1", "0.2")),
scale_x_continuous(limits = c(0, 0.3), breaks = c(0, 0.1, 0.2, 0.3), labels = c("0", "0.1", "0.2", "0.3")),
scale_x_continuous(limits = c(0, 0.5), breaks = c(0, 0.2, 0.4), labels = c("0", "0.2", "0.4"))
)
p_rmse_pointrange <- df %>%
filter(data_link_transformed == "Log") %>%
group_by(data_family_transformed, fit_family_transformed, shape) %>%
summarise(median_rmse_s = median(rmse_s),
min_rmse = quantile(rmse_s, probs = c(0.025)),
max_rmse = quantile(rmse_s, probs = c(0.975))) %>%
ggplot(
aes(y = data_family_transformed,
x = median_rmse_s,
xmin = min_rmse,
xmax = max_rmse,
color = fit_family_transformed)) +
geom_pointrange(size = 0.18, position = position_dodge(width = 0.8)) +
scale_color_scico_d(palette = "batlow", end = 0.8) +
facet_grid(~ shape, scales = "free_x") +
facetted_pos_scales(x = scales) +
xlab(TeX("$RMSE(beta_{xy})$")) +
ylab("Data likelihood") +
labs(color = "Fit likelihood") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
axis.ticks.y = element_blank())
## `summarise()` has grouped output by 'data_family_transformed',
## 'fit_family_transformed'. You can override using the `.groups` argument.
p_rmse_pointrange

ggsave(paste0(figure_path,"positive_rmse_showcase_pointrange.pdf"), width = 210, height = (297/4)*1.4, units = "mm", useDingbats = TRUE)
# RMSE softplus showcase
p_rmse_softplus_showcase <- df %>%
filter(data_link_transformed == "Softplus") %>%
ggplot(aes(y = fit_family_transformed, x = rmse_s)) +
coord_cartesian(xlim = c(0, 0.75)) +
scale_x_continuous(labels = c("0", "0.75"), breaks = c(0, 0.75)) +
geom_boxplot() +
stat_summary(fun="mean", shape = 4, size = 0.3) +
facet_grid(shape ~ data_family_transformed) +
xlab(TeX("$RMSE(beta_{xy})$")) +
ylab("Fit likelihood") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
axis.ticks.y = element_blank())
p_rmse_softplus_showcase
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).

ggsave(paste0(figure_path,"positive_rmse_softplus_showcase.pdf"), width = 210, height = (297/4)*1.5, units = "mm", useDingbats = TRUE)
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
# full absolute bias results
p_bias <- df %>%
ggplot(aes(y = fit_family_transformed, x = abs_bias)) +
coord_cartesian(xlim = c(0, 0.5)) +
scale_x_continuous(labels = c("0", "0.25", "0.5"), breaks = c(0, 0.25, 0.5)) +
geom_boxplot() +
stat_summary(fun="mean", shape = 4, size = 0.3) +
facet_grid(data_link_transformed + shape ~ data_family_transformed) +
xlab(TeX("$abs(bias(beta_{xy}))$")) +
ylab("Fit likelihood") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
axis.ticks.y = element_blank())
p_bias
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).

ggsave(paste0(figure_path, "positive_bias.pdf"), width = 210, height = (297/4) * 3, units = "mm", useDingbats = TRUE)
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
# absolute bias showcase
p_bias_showcase <- df %>%
filter(data_link_transformed == "Log") %>%
ggplot(aes(y = fit_family_transformed, x = abs_bias)) +
coord_cartesian(xlim = c(0, 0.25)) +
scale_x_continuous(labels = c("0", "0.25"), breaks = c(0, 0.25)) +
geom_boxplot() +
stat_summary(fun="mean", shape = 4, size = 0.3) +
facet_grid(shape ~ data_family_transformed) +
xlab(TeX("$abs(bias)(beta_{xy})$")) +
ylab("Fit likelihood") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
axis.ticks.y = element_blank())
p_bias_showcase
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).

ggsave(paste0(figure_path,"positive_bias_showcase.pdf"), width = 210, height = (297/4)*1.5, units = "mm", useDingbats = TRUE)
## Warning: Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
## Removed 7 rows containing missing values (`geom_segment()`).
ROC
calibration_df <- result_df %>%
# Calculate significance of the .
mutate(x_y_coef_cat = ifelse(x_y_coef == 0, 0, 1),
sig95 = pq_0.025 > 0 | pq_0.975 < 0,
sig90 = pq_0.05 > 0 | pq_0.95 < 0,
sig80 = pq_0.1 > 0 | pq_0.9 < 0,
sig50 = pq_0.25 > 0 | pq_0.75 < 0)
# We add the normal identity model to all fit links for easier visual comparison
calibration_df <- calibration_df %>%
mutate(
fit_family_transformed = case_when(
(fit_family_transformed == "Normal" & fit_link_transformed == "Identity") ~ "Normal Identity",
TRUE ~ fit_family_transformed)) %>%
mutate(
fit_link_transformed = case_when(
(fit_family_transformed == "Normal Identity") ~ "Log",
TRUE ~ fit_link_transformed))
data1 <- calibration_df %>%
filter(fit_family_transformed == "Normal Identity") %>%
mutate(fit_link_transformed = "Softplus")
calibration_df <- rbind(calibration_df, data1)
calibration_df = calibration_df %>%
mutate(fit_link_transformed = factor(fit_link_transformed, levels =
c("Log",
"Softplus")))
cols = c("formula", "fit_link_transformed", "fit_family_transformed", "shape", "data_link_transformed", "data_family_transformed")
TPR <- calibration_df %>%
filter(x_y_coef != 0) %>%
group_by(across(all_of(cols))) %>%
summarise(TPR_95 = mean(sig95),
TPR_90 = mean(sig90),
TPR_80 = mean(sig80),
TPR_50 = mean(sig50),
TPR_0 = 0,
TPR_1 = 1)
## `summarise()` has grouped output by 'formula', 'fit_link_transformed',
## 'fit_family_transformed', 'shape', 'data_link_transformed'. You can override
## using the `.groups` argument.
FPR <- calibration_df %>%
filter(x_y_coef == 0) %>%
group_by(across(all_of(cols))) %>%
summarise(FPR_95 = mean(sig95),
FPR_90 = mean(sig90),
FPR_80 = mean(sig80),
FPR_50 = mean(sig50),
FPR_0 = 0,
FPR_1 = 1)
## `summarise()` has grouped output by 'formula', 'fit_link_transformed',
## 'fit_family_transformed', 'shape', 'data_link_transformed'. You can override
## using the `.groups` argument.
df = merge(TPR, FPR)
auc <- vector(mode = "numeric", length = nrow(df))
for (i in seq_len(nrow(df))) {
auc[[i]] <- AUC(
x = c(df$FPR_0[[i]], df$FPR_95[[i]], df$FPR_90[[i]], df$FPR_80[[i]], df$FPR_50[[i]],df$FPR_1[[i]]),
y = c(df$TPR_0[[i]], df$TPR_95[[i]], df$TPR_90[[i]], df$TPR_80[[i]], df$TPR_50[[i]], df$TPR_1[[i]])
)
}
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
## Warning in regularize.values(x, y, ties, missing(ties), na.rm = na.rm):
## collapsing to unique 'x' values
df$auc <- auc
df <- df %>%
mutate(formula = sum_coding(
factor(formula,
levels = c("y ~ x + z1 + z2",
"y ~ x + z1",
"y ~ x + z1 + z2 + z3")))) %>%
mutate(data_link_transformed = sum_coding(
factor(data_link_transformed, levels = c("Log", "Softplus"))
)) %>%
mutate(fit_family_transformed = sum_coding(
factor(fit_family_transformed,
levels = c("Gamma",
"Weibull",
"Beta prime",
"trans. Normal",
"Fréchet",
"Gompertz",
"Normal",
"Normal Identity")))) %>%
mutate(shape = sum_coding(
factor(shape, levels = c("Thin Tail", "Heavy Tail", "Ramp")))
)
m1 = brm(auc ~ 1 +
fit_link_transformed + fit_family_transformed +
data_link_transformed + data_family_transformed +
shape + formula +
fit_link_transformed * data_link_transformed +
fit_family_transformed * data_family_transformed,
data = df, family = Beta(), cores = 4, warmup = 1000, iter = 3500,file = paste0(data_path, "m_auc.RDS"))
summary(m1)
## Family: beta
## Links: mu = logit; phi = identity
## Formula: auc ~ 1 + fit_link_transformed + fit_family_transformed + data_link_transformed + data_family_transformed + shape + formula + fit_link_transformed * data_link_transformed + fit_family_transformed * data_family_transformed
## Data: df (Number of observations: 1725)
## Draws: 4 chains, each with iter = 3500; warmup = 1000; thin = 1;
## total post-warmup draws = 10000
##
## Population-Level Effects:
## Estimate
## Intercept 2.02
## fit_link_transformedSoftplus -0.18
## fit_family_transformedWeibull 0.07
## fit_family_transformedBetaprime 0.04
## fit_family_transformedtrans.Normal 0.32
## fit_family_transformedFréchet -0.74
## fit_family_transformedGompertz -0.23
## fit_family_transformedNormal 0.01
## fit_family_transformedNormalIdentity 0.15
## data_link_transformedSoftplus -0.01
## data_family_transformedWeibull 0.14
## data_family_transformedFréchet 0.25
## data_family_transformedBetaprime 0.23
## data_family_transformedGompertz -0.09
## data_family_transformedtrans.Normal -0.32
## shapeHeavyTail 0.02
## shapeRamp 0.06
## formulay~xPz1 0.03
## formulay~xPz1Pz2Pz3 -0.24
## fit_link_transformedSoftplus:data_link_transformedSoftplus 0.20
## fit_family_transformedWeibull:data_family_transformedWeibull 0.27
## fit_family_transformedBetaprime:data_family_transformedWeibull -0.26
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 0.00
## fit_family_transformedFréchet:data_family_transformedWeibull -0.45
## fit_family_transformedGompertz:data_family_transformedWeibull 0.25
## fit_family_transformedNormal:data_family_transformedWeibull 0.04
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 0.10
## fit_family_transformedWeibull:data_family_transformedFréchet -0.85
## fit_family_transformedBetaprime:data_family_transformedFréchet 0.75
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 0.10
## fit_family_transformedFréchet:data_family_transformedFréchet 1.50
## fit_family_transformedGompertz:data_family_transformedFréchet -0.66
## fit_family_transformedNormal:data_family_transformedFréchet -0.17
## fit_family_transformedNormalIdentity:data_family_transformedFréchet -0.24
## fit_family_transformedWeibull:data_family_transformedBetaprime -0.34
## fit_family_transformedBetaprime:data_family_transformedBetaprime 0.45
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 0.19
## fit_family_transformedFréchet:data_family_transformedBetaprime 0.46
## fit_family_transformedGompertz:data_family_transformedBetaprime -0.37
## fit_family_transformedNormal:data_family_transformedBetaprime -0.03
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime -0.08
## fit_family_transformedWeibull:data_family_transformedGompertz 0.48
## fit_family_transformedBetaprime:data_family_transformedGompertz -0.83
## fit_family_transformedtrans.Normal:data_family_transformedGompertz -0.11
## fit_family_transformedFréchet:data_family_transformedGompertz -0.75
## fit_family_transformedGompertz:data_family_transformedGompertz 0.71
## fit_family_transformedNormal:data_family_transformedGompertz 0.30
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 0.25
## fit_family_transformedWeibull:data_family_transformedtrans.Normal -0.18
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 0.17
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 0.01
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 0.37
## fit_family_transformedGompertz:data_family_transformedtrans.Normal -0.06
## fit_family_transformedNormal:data_family_transformedtrans.Normal -0.00
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal -0.06
## Est.Error
## Intercept 0.03
## fit_link_transformedSoftplus 0.03
## fit_family_transformedWeibull 0.08
## fit_family_transformedBetaprime 0.08
## fit_family_transformedtrans.Normal 0.09
## fit_family_transformedFréchet 0.07
## fit_family_transformedGompertz 0.08
## fit_family_transformedNormal 0.08
## fit_family_transformedNormalIdentity 0.08
## data_link_transformedSoftplus 0.02
## data_family_transformedWeibull 0.05
## data_family_transformedFréchet 0.05
## data_family_transformedBetaprime 0.05
## data_family_transformedGompertz 0.04
## data_family_transformedtrans.Normal 0.04
## shapeHeavyTail 0.02
## shapeRamp 0.02
## formulay~xPz1 0.02
## formulay~xPz1Pz2Pz3 0.02
## fit_link_transformedSoftplus:data_link_transformedSoftplus 0.03
## fit_family_transformedWeibull:data_family_transformedWeibull 0.13
## fit_family_transformedBetaprime:data_family_transformedWeibull 0.11
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 0.13
## fit_family_transformedFréchet:data_family_transformedWeibull 0.10
## fit_family_transformedGompertz:data_family_transformedWeibull 0.11
## fit_family_transformedNormal:data_family_transformedWeibull 0.12
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 0.12
## fit_family_transformedWeibull:data_family_transformedFréchet 0.11
## fit_family_transformedBetaprime:data_family_transformedFréchet 0.14
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 0.14
## fit_family_transformedFréchet:data_family_transformedFréchet 0.13
## fit_family_transformedGompertz:data_family_transformedFréchet 0.11
## fit_family_transformedNormal:data_family_transformedFréchet 0.12
## fit_family_transformedNormalIdentity:data_family_transformedFréchet 0.12
## fit_family_transformedWeibull:data_family_transformedBetaprime 0.12
## fit_family_transformedBetaprime:data_family_transformedBetaprime 0.13
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 0.14
## fit_family_transformedFréchet:data_family_transformedBetaprime 0.11
## fit_family_transformedGompertz:data_family_transformedBetaprime 0.11
## fit_family_transformedNormal:data_family_transformedBetaprime 0.12
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime 0.12
## fit_family_transformedWeibull:data_family_transformedGompertz 0.12
## fit_family_transformedBetaprime:data_family_transformedGompertz 0.11
## fit_family_transformedtrans.Normal:data_family_transformedGompertz 0.12
## fit_family_transformedFréchet:data_family_transformedGompertz 0.09
## fit_family_transformedGompertz:data_family_transformedGompertz 0.12
## fit_family_transformedNormal:data_family_transformedGompertz 0.12
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 0.12
## fit_family_transformedWeibull:data_family_transformedtrans.Normal 0.11
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 0.11
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 0.12
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 0.10
## fit_family_transformedGompertz:data_family_transformedtrans.Normal 0.10
## fit_family_transformedNormal:data_family_transformedtrans.Normal 0.11
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal 0.11
## l-95% CI
## Intercept 1.95
## fit_link_transformedSoftplus -0.23
## fit_family_transformedWeibull -0.09
## fit_family_transformedBetaprime -0.12
## fit_family_transformedtrans.Normal 0.15
## fit_family_transformedFréchet -0.88
## fit_family_transformedGompertz -0.38
## fit_family_transformedNormal -0.15
## fit_family_transformedNormalIdentity -0.01
## data_link_transformedSoftplus -0.04
## data_family_transformedWeibull 0.05
## data_family_transformedFréchet 0.16
## data_family_transformedBetaprime 0.14
## data_family_transformedGompertz -0.18
## data_family_transformedtrans.Normal -0.41
## shapeHeavyTail -0.02
## shapeRamp 0.03
## formulay~xPz1 -0.01
## formulay~xPz1Pz2Pz3 -0.27
## fit_link_transformedSoftplus:data_link_transformedSoftplus 0.15
## fit_family_transformedWeibull:data_family_transformedWeibull 0.03
## fit_family_transformedBetaprime:data_family_transformedWeibull -0.48
## fit_family_transformedtrans.Normal:data_family_transformedWeibull -0.25
## fit_family_transformedFréchet:data_family_transformedWeibull -0.64
## fit_family_transformedGompertz:data_family_transformedWeibull 0.03
## fit_family_transformedNormal:data_family_transformedWeibull -0.20
## fit_family_transformedNormalIdentity:data_family_transformedWeibull -0.14
## fit_family_transformedWeibull:data_family_transformedFréchet -1.07
## fit_family_transformedBetaprime:data_family_transformedFréchet 0.49
## fit_family_transformedtrans.Normal:data_family_transformedFréchet -0.16
## fit_family_transformedFréchet:data_family_transformedFréchet 1.26
## fit_family_transformedGompertz:data_family_transformedFréchet -0.86
## fit_family_transformedNormal:data_family_transformedFréchet -0.41
## fit_family_transformedNormalIdentity:data_family_transformedFréchet -0.48
## fit_family_transformedWeibull:data_family_transformedBetaprime -0.56
## fit_family_transformedBetaprime:data_family_transformedBetaprime 0.19
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime -0.08
## fit_family_transformedFréchet:data_family_transformedBetaprime 0.25
## fit_family_transformedGompertz:data_family_transformedBetaprime -0.58
## fit_family_transformedNormal:data_family_transformedBetaprime -0.27
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime -0.32
## fit_family_transformedWeibull:data_family_transformedGompertz 0.24
## fit_family_transformedBetaprime:data_family_transformedGompertz -1.04
## fit_family_transformedtrans.Normal:data_family_transformedGompertz -0.35
## fit_family_transformedFréchet:data_family_transformedGompertz -0.94
## fit_family_transformedGompertz:data_family_transformedGompertz 0.48
## fit_family_transformedNormal:data_family_transformedGompertz 0.06
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 0.01
## fit_family_transformedWeibull:data_family_transformedtrans.Normal -0.39
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal -0.05
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal -0.23
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 0.17
## fit_family_transformedGompertz:data_family_transformedtrans.Normal -0.27
## fit_family_transformedNormal:data_family_transformedtrans.Normal -0.22
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal -0.28
## u-95% CI
## Intercept 2.08
## fit_link_transformedSoftplus -0.13
## fit_family_transformedWeibull 0.24
## fit_family_transformedBetaprime 0.21
## fit_family_transformedtrans.Normal 0.50
## fit_family_transformedFréchet -0.60
## fit_family_transformedGompertz -0.07
## fit_family_transformedNormal 0.18
## fit_family_transformedNormalIdentity 0.32
## data_link_transformedSoftplus 0.03
## data_family_transformedWeibull 0.23
## data_family_transformedFréchet 0.34
## data_family_transformedBetaprime 0.32
## data_family_transformedGompertz -0.00
## data_family_transformedtrans.Normal -0.24
## shapeHeavyTail 0.05
## shapeRamp 0.10
## formulay~xPz1 0.06
## formulay~xPz1Pz2Pz3 -0.20
## fit_link_transformedSoftplus:data_link_transformedSoftplus 0.25
## fit_family_transformedWeibull:data_family_transformedWeibull 0.53
## fit_family_transformedBetaprime:data_family_transformedWeibull -0.03
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 0.26
## fit_family_transformedFréchet:data_family_transformedWeibull -0.26
## fit_family_transformedGompertz:data_family_transformedWeibull 0.48
## fit_family_transformedNormal:data_family_transformedWeibull 0.28
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 0.34
## fit_family_transformedWeibull:data_family_transformedFréchet -0.63
## fit_family_transformedBetaprime:data_family_transformedFréchet 1.03
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 0.36
## fit_family_transformedFréchet:data_family_transformedFréchet 1.76
## fit_family_transformedGompertz:data_family_transformedFréchet -0.45
## fit_family_transformedNormal:data_family_transformedFréchet 0.07
## fit_family_transformedNormalIdentity:data_family_transformedFréchet -0.00
## fit_family_transformedWeibull:data_family_transformedBetaprime -0.11
## fit_family_transformedBetaprime:data_family_transformedBetaprime 0.71
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 0.46
## fit_family_transformedFréchet:data_family_transformedBetaprime 0.67
## fit_family_transformedGompertz:data_family_transformedBetaprime -0.15
## fit_family_transformedNormal:data_family_transformedBetaprime 0.21
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime 0.16
## fit_family_transformedWeibull:data_family_transformedGompertz 0.72
## fit_family_transformedBetaprime:data_family_transformedGompertz -0.61
## fit_family_transformedtrans.Normal:data_family_transformedGompertz 0.13
## fit_family_transformedFréchet:data_family_transformedGompertz -0.57
## fit_family_transformedGompertz:data_family_transformedGompertz 0.95
## fit_family_transformedNormal:data_family_transformedGompertz 0.54
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 0.49
## fit_family_transformedWeibull:data_family_transformedtrans.Normal 0.04
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 0.39
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 0.24
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 0.56
## fit_family_transformedGompertz:data_family_transformedtrans.Normal 0.14
## fit_family_transformedNormal:data_family_transformedtrans.Normal 0.22
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal 0.17
## Rhat
## Intercept 1.00
## fit_link_transformedSoftplus 1.00
## fit_family_transformedWeibull 1.00
## fit_family_transformedBetaprime 1.00
## fit_family_transformedtrans.Normal 1.00
## fit_family_transformedFréchet 1.00
## fit_family_transformedGompertz 1.00
## fit_family_transformedNormal 1.00
## fit_family_transformedNormalIdentity 1.00
## data_link_transformedSoftplus 1.00
## data_family_transformedWeibull 1.00
## data_family_transformedFréchet 1.00
## data_family_transformedBetaprime 1.00
## data_family_transformedGompertz 1.00
## data_family_transformedtrans.Normal 1.00
## shapeHeavyTail 1.00
## shapeRamp 1.00
## formulay~xPz1 1.00
## formulay~xPz1Pz2Pz3 1.00
## fit_link_transformedSoftplus:data_link_transformedSoftplus 1.00
## fit_family_transformedWeibull:data_family_transformedWeibull 1.00
## fit_family_transformedBetaprime:data_family_transformedWeibull 1.00
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 1.00
## fit_family_transformedFréchet:data_family_transformedWeibull 1.00
## fit_family_transformedGompertz:data_family_transformedWeibull 1.00
## fit_family_transformedNormal:data_family_transformedWeibull 1.00
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 1.00
## fit_family_transformedWeibull:data_family_transformedFréchet 1.00
## fit_family_transformedBetaprime:data_family_transformedFréchet 1.00
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 1.00
## fit_family_transformedFréchet:data_family_transformedFréchet 1.00
## fit_family_transformedGompertz:data_family_transformedFréchet 1.00
## fit_family_transformedNormal:data_family_transformedFréchet 1.00
## fit_family_transformedNormalIdentity:data_family_transformedFréchet 1.00
## fit_family_transformedWeibull:data_family_transformedBetaprime 1.00
## fit_family_transformedBetaprime:data_family_transformedBetaprime 1.00
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 1.00
## fit_family_transformedFréchet:data_family_transformedBetaprime 1.00
## fit_family_transformedGompertz:data_family_transformedBetaprime 1.00
## fit_family_transformedNormal:data_family_transformedBetaprime 1.00
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime 1.00
## fit_family_transformedWeibull:data_family_transformedGompertz 1.00
## fit_family_transformedBetaprime:data_family_transformedGompertz 1.00
## fit_family_transformedtrans.Normal:data_family_transformedGompertz 1.00
## fit_family_transformedFréchet:data_family_transformedGompertz 1.00
## fit_family_transformedGompertz:data_family_transformedGompertz 1.00
## fit_family_transformedNormal:data_family_transformedGompertz 1.00
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 1.00
## fit_family_transformedWeibull:data_family_transformedtrans.Normal 1.00
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 1.00
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 1.00
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 1.00
## fit_family_transformedGompertz:data_family_transformedtrans.Normal 1.00
## fit_family_transformedNormal:data_family_transformedtrans.Normal 1.00
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal 1.00
## Bulk_ESS
## Intercept 7384
## fit_link_transformedSoftplus 21790
## fit_family_transformedWeibull 4657
## fit_family_transformedBetaprime 4906
## fit_family_transformedtrans.Normal 4441
## fit_family_transformedFréchet 4603
## fit_family_transformedGompertz 4574
## fit_family_transformedNormal 4065
## fit_family_transformedNormalIdentity 3941
## data_link_transformedSoftplus 12494
## data_family_transformedWeibull 8466
## data_family_transformedFréchet 8607
## data_family_transformedBetaprime 8579
## data_family_transformedGompertz 8180
## data_family_transformedtrans.Normal 8195
## shapeHeavyTail 14724
## shapeRamp 14513
## formulay~xPz1 14542
## formulay~xPz1Pz2Pz3 14188
## fit_link_transformedSoftplus:data_link_transformedSoftplus 12532
## fit_family_transformedWeibull:data_family_transformedWeibull 7124
## fit_family_transformedBetaprime:data_family_transformedWeibull 6807
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 6346
## fit_family_transformedFréchet:data_family_transformedWeibull 6061
## fit_family_transformedGompertz:data_family_transformedWeibull 6804
## fit_family_transformedNormal:data_family_transformedWeibull 6322
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 6180
## fit_family_transformedWeibull:data_family_transformedFréchet 6133
## fit_family_transformedBetaprime:data_family_transformedFréchet 7679
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 6557
## fit_family_transformedFréchet:data_family_transformedFréchet 8506
## fit_family_transformedGompertz:data_family_transformedFréchet 6195
## fit_family_transformedNormal:data_family_transformedFréchet 5794
## fit_family_transformedNormalIdentity:data_family_transformedFréchet 5763
## fit_family_transformedWeibull:data_family_transformedBetaprime 6816
## fit_family_transformedBetaprime:data_family_transformedBetaprime 7490
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 6584
## fit_family_transformedFréchet:data_family_transformedBetaprime 6965
## fit_family_transformedGompertz:data_family_transformedBetaprime 6173
## fit_family_transformedNormal:data_family_transformedBetaprime 5845
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime 5800
## fit_family_transformedWeibull:data_family_transformedGompertz 6584
## fit_family_transformedBetaprime:data_family_transformedGompertz 6180
## fit_family_transformedtrans.Normal:data_family_transformedGompertz 5728
## fit_family_transformedFréchet:data_family_transformedGompertz 5659
## fit_family_transformedGompertz:data_family_transformedGompertz 6682
## fit_family_transformedNormal:data_family_transformedGompertz 5793
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 5753
## fit_family_transformedWeibull:data_family_transformedtrans.Normal 6194
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 6704
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 5629
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 6547
## fit_family_transformedGompertz:data_family_transformedtrans.Normal 6221
## fit_family_transformedNormal:data_family_transformedtrans.Normal 5538
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal 5426
## Tail_ESS
## Intercept 7596
## fit_link_transformedSoftplus 6743
## fit_family_transformedWeibull 6266
## fit_family_transformedBetaprime 6607
## fit_family_transformedtrans.Normal 5942
## fit_family_transformedFréchet 6139
## fit_family_transformedGompertz 6723
## fit_family_transformedNormal 5858
## fit_family_transformedNormalIdentity 5473
## data_link_transformedSoftplus 8009
## data_family_transformedWeibull 8234
## data_family_transformedFréchet 8317
## data_family_transformedBetaprime 8610
## data_family_transformedGompertz 7766
## data_family_transformedtrans.Normal 8376
## shapeHeavyTail 8467
## shapeRamp 8341
## formulay~xPz1 8046
## formulay~xPz1Pz2Pz3 7669
## fit_link_transformedSoftplus:data_link_transformedSoftplus 7717
## fit_family_transformedWeibull:data_family_transformedWeibull 7559
## fit_family_transformedBetaprime:data_family_transformedWeibull 7822
## fit_family_transformedtrans.Normal:data_family_transformedWeibull 7340
## fit_family_transformedFréchet:data_family_transformedWeibull 6523
## fit_family_transformedGompertz:data_family_transformedWeibull 7973
## fit_family_transformedNormal:data_family_transformedWeibull 7334
## fit_family_transformedNormalIdentity:data_family_transformedWeibull 7683
## fit_family_transformedWeibull:data_family_transformedFréchet 7432
## fit_family_transformedBetaprime:data_family_transformedFréchet 7719
## fit_family_transformedtrans.Normal:data_family_transformedFréchet 7617
## fit_family_transformedFréchet:data_family_transformedFréchet 7418
## fit_family_transformedGompertz:data_family_transformedFréchet 7562
## fit_family_transformedNormal:data_family_transformedFréchet 6664
## fit_family_transformedNormalIdentity:data_family_transformedFréchet 6805
## fit_family_transformedWeibull:data_family_transformedBetaprime 7940
## fit_family_transformedBetaprime:data_family_transformedBetaprime 7748
## fit_family_transformedtrans.Normal:data_family_transformedBetaprime 7344
## fit_family_transformedFréchet:data_family_transformedBetaprime 6998
## fit_family_transformedGompertz:data_family_transformedBetaprime 7533
## fit_family_transformedNormal:data_family_transformedBetaprime 7010
## fit_family_transformedNormalIdentity:data_family_transformedBetaprime 6974
## fit_family_transformedWeibull:data_family_transformedGompertz 7276
## fit_family_transformedBetaprime:data_family_transformedGompertz 7153
## fit_family_transformedtrans.Normal:data_family_transformedGompertz 7177
## fit_family_transformedFréchet:data_family_transformedGompertz 7340
## fit_family_transformedGompertz:data_family_transformedGompertz 7700
## fit_family_transformedNormal:data_family_transformedGompertz 6951
## fit_family_transformedNormalIdentity:data_family_transformedGompertz 7054
## fit_family_transformedWeibull:data_family_transformedtrans.Normal 7304
## fit_family_transformedBetaprime:data_family_transformedtrans.Normal 8020
## fit_family_transformedtrans.Normal:data_family_transformedtrans.Normal 7629
## fit_family_transformedFréchet:data_family_transformedtrans.Normal 7461
## fit_family_transformedGompertz:data_family_transformedtrans.Normal 7605
## fit_family_transformedNormal:data_family_transformedtrans.Normal 6990
## fit_family_transformedNormalIdentity:data_family_transformedtrans.Normal 6783
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## phi 28.61 0.99 26.72 30.59 1.00 16296 8310
##
## Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
p1 = conditional_effects(m1, effects = "data_link_transformed:fit_link_transformed",
conditions = data.frame(
shape = NA,
data_family_transformed = NA,
fit_family_transformed = NA,
formula = NA)
)
plot(p1, plot = FALSE)[[1]] + scale_color_scico_d(palette = "batlow", end = 0.8) +
xlab("Data Link") + ylab("AUC")+
theme_bw(base_size = 12) +
labs(color = "Fit Link", fill = "Fit Link")

ggsave(paste0(figure_path,"positive_auc_cond_link_link.pdf"), width = 210, height = (297/4)/1.5, units = "mm", useDingbats = TRUE)
p1$`data_link_transformed:fit_link_transformed` %>%
ggplot(
aes(y = estimate__,
x = effect1__,
ymin = lower__,
ymax = upper__,
color = effect2__)) +
geom_pointrange(size = 0.4, linewidth = 0.8, position = position_dodge(width = 0.2)) +
scale_color_scico_d(palette = "batlow", end = 0.8) +
xlab("Data Link") +
ylab("AUC ROC") +
labs(color = "Fit Link") +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 9),
strip.text.y = element_text(size = 9),
axis.ticks.y = element_blank())

ggsave(paste0(figure_path,"positive_auc_cond_link_link.pdf"), width = 210, height = (297/4)*0.6, units = "mm", useDingbats = TRUE)
p2 = conditional_effects(m1, effects = "data_family_transformed:fit_family_transformed",
conditions = data.frame(
shape = NA,
data_link_transformed = NA,
fit_link_transformed = NA,
formula = NA))
plot(p2, plot = FALSE)[[1]] + scale_color_scico_d(palette = "batlow", end = 0.8) +
xlab("Data Family") + ylab("AUC")+
theme_bw(base_size = 12) +
labs(color = "Fit Family", fill = "Fit Family")

ggsave(paste0(figure_path,"positive_auc_cond_fam_fam.pdf"), width = 210, height = (297/4)*1.2, units = "mm", useDingbats = TRUE)
p3 = conditional_effects(m1, effects = "fit_link_transformed",
conditions = data.frame(
shape = NA,
data_family_transformed = NA,
fit_family_transformed = NA,
data_link_transformed = NA,
formula = NA)
)
plot(p3, plot = FALSE)[[1]] + scale_color_scico_d(palette = "batlow", end = 0.8) +
xlab("Fit Link") + ylab("AUC")+
theme_bw(base_size = 12)

ggsave(paste0(figure_path,"positive_auc_cond_link.pdf"), width = 210/2, height = (297/4)*1.2, units = "mm", useDingbats = TRUE)
p4 = conditional_effects(m1, effects = "fit_family_transformed",
conditions = data.frame(
shape = NA,
data_family_transformed = NA,
fit_link_transformed = NA,
data_link_transformed = NA,
formula = NA)
)
plot(p4, plot = FALSE)[[1]] + scale_color_scico_d(palette = "batlow", end = 0.8) +
xlab("Fit Family") + ylab("AUC")+
theme_bw(base_size = 12)

ggsave(paste0(figure_path,"positive_auc_cond_fam.pdf"), width = 210, height = (297/4)*1.2, units = "mm", useDingbats = TRUE)
# Function to prepare a TPR/FPR DF averaging over different columns
roc_df <- function(cols){
TPR <- calibration_df %>%
filter(x_y_coef != 0) %>%
group_by(across(all_of(cols))) %>%
summarise(TPR_95 = mean(sig95),
TPR_90 = mean(sig90),
TPR_80 = mean(sig80),
TPR_50 = mean(sig50),
TPR_0 = 0,
TPR_1 = 1)
FPR <- calibration_df %>%
filter(x_y_coef == 0) %>%
group_by(across(all_of(cols))) %>%
summarise(FPR_95 = mean(sig95),
FPR_90 = mean(sig90),
FPR_80 = mean(sig80),
FPR_50 = mean(sig50),
FPR_0 = 0,
FPR_1 = 1)
df = merge(TPR, FPR)
df = df %>%
pivot_longer(
cols = all_of(colnames(df)[! colnames(df) %in% cols]),
names_pattern = "TPR_(.*)$|FPR_(.*)$",
names_to = c("TPR", "FPR")) %>%
mutate(TPR = as.numeric(TPR),
FPR = as.numeric(FPR)) %>%
mutate(CI = if_else(is.na(TPR), FPR, TPR),
TPR = if_else(is.na(TPR), NA, value),
FPR = if_else(is.na(FPR), NA, value))
part1 <- filter(df, is.na(TPR))
part1$TPR <- NULL
part1$value <- NULL
part2 <- filter(df, is.na(FPR))
part2$FPR <- NULL
part2$value <- NULL
df = merge(part1, part2)
return(df)
}
# Full ROC Plot
df <- roc_df(cols = c("formula", "fit_link_transformed", "fit_family_transformed", "shape", "data_link_transformed", "data_family_transformed"))
## `summarise()` has grouped output by 'formula', 'fit_link_transformed',
## 'fit_family_transformed', 'shape', 'data_link_transformed'. You can override
## using the `.groups` argument.
## `summarise()` has grouped output by 'formula', 'fit_link_transformed',
## 'fit_family_transformed', 'shape', 'data_link_transformed'. You can override
## using the `.groups` argument.
roc <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1),labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed), data = filter(df, CI == "95"), size = 2) +
facet_grid(shape + data_link_transformed + data_family_transformed ~ formula + fit_link_transformed) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc

ggsave(paste0(figure_path, "positive_roc.pdf"), width = 210*2, height = (297/4)*16, units = "mm", useDingbats = TRUE)
# Compare the fit configurations and average over all DGPS
df <- roc_df(cols = c("formula", "fit_link_transformed", "fit_family_transformed"))
## `summarise()` has grouped output by 'formula', 'fit_link_transformed'. You can
## override using the `.groups` argument.
## `summarise()` has grouped output by 'formula', 'fit_link_transformed'. You can
## override using the `.groups` argument.
roc_by_fit <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed), data = filter(df, CI == "95"), size = 2) +
facet_grid(fit_link_transformed ~ formula) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_fit

ggsave(paste0(figure_path, "positive_roc_by_fit.pdf"), width = 210, height = (297/4) * 1.5, units = "mm", useDingbats = TRUE)
# Compare the DGPs and average over all fit configurations
df <- roc_df(cols = c("shape", "data_link_transformed", "data_family_transformed", "fit_family_transformed"))
## `summarise()` has grouped output by 'shape', 'data_link_transformed',
## 'data_family_transformed'. You can override using the `.groups` argument.
## `summarise()` has grouped output by 'shape', 'data_link_transformed',
## 'data_family_transformed'. You can override using the `.groups` argument.
roc_by_dgp <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed), data = filter(df, CI == "95"), size = 2) +
facet_grid(shape + data_link_transformed ~ data_family_transformed) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_dgp

ggsave(paste0(figure_path, "positive_roc_by_dgp.pdf"), width = 210, height = (297/4) * 3, units = "mm", useDingbats = TRUE)
# Compare the DGP shape and fit configurations
df <- roc_df(cols = c("fit_link_transformed", "fit_family_transformed", "shape"))
## `summarise()` has grouped output by 'fit_link_transformed',
## 'fit_family_transformed'. You can override using the `.groups` argument.
## `summarise()` has grouped output by 'fit_link_transformed',
## 'fit_family_transformed'. You can override using the `.groups` argument.
roc_by_shape_and_fit <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed),
data = filter(df, CI == "95"), size = 2) +
#geom_point(aes(color = fit_family_transformed), alpha = 0.5) +
facet_grid(fit_link_transformed ~ shape) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_shape_and_fit

ggsave(paste0(figure_path, "positive_roc_shape_and_fit.pdf"), width = 210, height = (297/4) * 1.2, units = "mm", useDingbats = TRUE)
# Compare the DGP Likelihood and fit configurations
df <- roc_df(cols = c("fit_link_transformed", "fit_family_transformed", "data_family_transformed"))
## `summarise()` has grouped output by 'fit_link_transformed',
## 'fit_family_transformed'. You can override using the `.groups` argument.
## `summarise()` has grouped output by 'fit_link_transformed',
## 'fit_family_transformed'. You can override using the `.groups` argument.
roc_by_shape_and_fit <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed),
data = filter(df, CI == "95"), size = 2) +
#geom_point(aes(color = fit_family_transformed), alpha = 0.5) +
facet_grid(fit_link_transformed ~ data_family_transformed) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_shape_and_fit

ggsave(paste0(figure_path, "positive_roc_dgp_lik_and_fit.pdf"), width = 210, height = (297/4) * 1.2, units = "mm", useDingbats = TRUE)
# Compare the links
df <- roc_df(cols = c("fit_link_transformed", "data_link_transformed"))
## `summarise()` has grouped output by 'fit_link_transformed'. You can override
## using the `.groups` argument.
## `summarise()` has grouped output by 'fit_link_transformed'. You can override
## using the `.groups` argument.
roc_by_shape_and_fit <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_link_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_link_transformed),
data = filter(df, CI == "95"), size = 2) +
#geom_point(aes(color = fit_family_transformed), alpha = 0.5) +
facet_grid(. ~ data_link_transformed) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_shape_and_fit

ggsave(paste0(figure_path, "positive_roc_links.png"), width = 210*2, height = (297/4) * 1.2*2, units = "mm")#, useDingbats = TRUE)
# Compare the likelihoods
df <- roc_df(cols = c("fit_family_transformed", "data_family_transformed"))
## `summarise()` has grouped output by 'fit_family_transformed'. You can override
## using the `.groups` argument.
## `summarise()` has grouped output by 'fit_family_transformed'. You can override
## using the `.groups` argument.
roc_by_shape_and_fit <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_family_transformed), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_family_transformed),
data = filter(df, CI == "95"), size = 2) +
#geom_point(aes(color = fit_family_transformed), alpha = 0.5) +
facet_grid(. ~ data_family_transformed) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=c(0, 1, 8, 4, 5, 17, 6, 7)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_shape_and_fit

ggsave(paste0(figure_path, "positive_roc_likelihoods.png"), width = 210*2, height = (297/4) * 1.2*2, units = "mm")#, useDingbats = TRUE)
# Compare LLs
df <- roc_df(cols = c("fit_family_transformed", "fit_link_transformed", "data_family_transformed", "data_link_transformed")) %>%
mutate(data_LL = paste0(data_family_transformed, data_link_transformed),
fit_LL = paste0(fit_family_transformed, fit_link_transformed))
## `summarise()` has grouped output by 'fit_family_transformed',
## 'fit_link_transformed', 'data_family_transformed'. You can override using the
## `.groups` argument.
## `summarise()` has grouped output by 'fit_family_transformed',
## 'fit_link_transformed', 'data_family_transformed'. You can override using the
## `.groups` argument.
roc_by_LL <- df %>%
ggplot(aes(FPR, TPR)) +
coord_cartesian(xlim = c(0,1), ylim = c(0,1)) +
scale_x_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
scale_y_continuous(breaks = c(0, 0.5, 1), labels = c("0", "0.5", "1")) +
geom_line(aes(color = fit_LL), linetype = "longdash", linewidth = 1) +
geom_abline(color = "grey", linetype = "dotted", slope = 1) +
geom_point(aes(FPR, TPR, shape = fit_LL),
data = filter(df, CI == "95"), size = 2) +
#geom_point(aes(color = fit_family_transformed), alpha = 0.5) +
facet_wrap(. ~ data_LL) +
scale_colour_scico_d(palette = "batlow", end = 0.8) +
scale_shape_manual(values=seq(1,16)) +
theme_bw(base_size = 12) +
theme(strip.text.x = element_text(size = 10),
strip.text.y = element_text(size = 10),
legend.position = "bottom") +
labs(color = "Fit family", shape = "Fit family") +
xlab("False positive rate") + ylab("True positive rate")
roc_by_LL

ggsave(paste0(figure_path, "positive_roc_LL.png"), width = 210*2, height = (297/4) * 1.2*3, units = "mm")#, useDingbats = TRUE)